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Abstract. The Large Area Telescope (LAT) instrument on the Gamma Ray Large Area Space 
Telescope (GLAST) has been designed to detect high-energy gamma rays and 'determine their 
direction of incidence and energy- We propose a reconstruction algorithm based on recent advances 
in statistical methodology. This method, alternative to the standard event analysis inherited from 
high energy collider physics experiments, incorporates more accurately the physical processes 
occurring in the detector, and makes full use of the statistical information available. It could thus 
provide a better estimate of the direction and energy of the primary photon. 


INTRODUCTION 

Gamma rays are produced by some of the highest energy events in the universe, and 
their study is critical to understanding the source locations and production mechanisms 
of ultra-relativistic cosmic particles [1, and refs therein]. 

The LAT instrument on the Gamma Ray Large Area Space Telescope (GLAST) 
satellite is scheduled for launch in late 2006. Its goal is to observe the gamma ray sky in 
the energy range 20 MeV to 300 GeV. 

At these energies, the photon interacts with matter primarily by annihilating into an 
electron positron pair. As a consequence, its parameters of interest (direction and energy) 
can be estimated only through the subsequent interaction of the two charged tracks 
with the active material of the detector. In the next section we describe in more detail 
this detection process and emphasize the challenges that it raises. We then outline our 
analysis methodology and show its application to a simulation of a simplified detector. 


THE LARGE AREA TELESCOPE INSTRUMENT 

The Large Area Telescope (LAT) instrument (Figure 1) is designed to determine the 
energy and direction of gamma ray photons that arrive at the detector. The detector is 
a stack of conversion and detection layers. The gamma ray photon is converted into 
an electron/positron pair in one of the conversion layers, and these particles cascade 
through the detector. The electron and positron are themselves subject to secondary 
scattering, and also produce further electrons and photons. All these particles trigger 
the detector strips, which read out the position where a particle crossed the detector 
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FIGURE 1. Exploded view of the LAT showing an incident gamma ray entering one tower. 


layer. The particles finally come to rest in the calorimeter, which measures the remaining 
energy. The anticoincidence detector eliminates event triggers due to charged cosmic 
rays (incident charged particles trigger the anticoincidence detector; incident photons 
do not.) Figure 1 shows a schematic of the LAT. The main output from the detector 
is the set of positions of the microstrip detectors that were triggered, together with the 
response of the calorimeter. It is from this set of positions and the calorimeter response 
that the energy and direction of the photon must be determined. 

Figure 2 shows a simulation of a 200 MeV photon arriving at the detector, obtained 
with the current software written by the GLAST collaboration [2]. The diagonal line 
extending to the top of the figure (blue) is the incident gamma ray, which and after 
conversion, splits into the two tracks followed by the electron and positron. It can be 
seen that the tracks are not straight, indicating significant secondary scattering. The 
small oblique squares (orange) indicate the microstrip detectors that fired during the 
event. The small cubes towards the bottom of the figure (blue and red) are the response 
of the calorimeter. The large rectangles on the right hand side and front (green) are tiles 
of the anticoincidence detector firing due to the incidence of a secondary electron. The 
tracks of some of the secondary particles are also shown (black lines). It is clear that 
most of the firings are due to the passage of the electron/positron pair, but that some 
are caused by hits from secondary electrons. These firings, together with the significant 
secondary scattering, represent a challenge for the correct reconstruction of the photon 
incident direction and energy. 


THE PHYSICS OF THE DETECTOR 

The primary physics process is pair production, where a gamma ray is converted into 
an electron/positron pair. The tungsten foils in the detector are dimensioned such that 
the conversion probability in a given foil is approximately 2% (and depends only very 



FIGURE 2. Simulation of the interaction of a 200MeV photon with the LAT. 


weakly on energy). This represents a design compromise - thicker foils would give 
greater conversion probability, but would also increase the effects of multiple scattering 
and other secondary physical processes which obscure the primary event. 

In the conversion event the energy of the incident photon is divided stochastically 
between the electron and positron. Figure 3 (left) shows the distribution of the energy 
split as a function of photon energy [3], Note that there is a reasonable probability that 
one of the secondary particles will be of very low energy, and thus very susceptible to 
being scattered or annihilated. The divergence angles of the electron and positron depend 
stochastically on the energies, as 
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where u is distributed as shown in figure 3 (right). 

The other physics processes that must be taken into account in the event analysis are: 

• multiple Coulomb scattering is the cumulative effect of multiple scatters of a 
charged particle traversing a medium. The name derives from the fact that the main 
component is the deflection due to Coulomb scattering on the nuclei. If we neglect 
the tails (corresponding to large scatter), the distribution of the deflection angle 
projected on a plane is approximately Gaussian, with a width given by: 

. 13.6MeV r— \ x . 

9 ° = "~p C p - ^W+o-omos-] 




0.6, 



pholOfl «n*<gy (0«V) 


FIGURE 3. Distribution of the energy split between the electron and positron (left). Distribution of the 
auxiliary quantity in determining the electron/positron angles (right). 



where p, fic and z are the momentum, velocity, and charge number of the incident 
particle, and x/X Q is the thickness of the scattering medium in radiation lengths [3], 

• positron annihilation is the interaction of a positron with an electron of the 
medium, resulting in the annihilation of the positron and electron, and the creation 
of two photons. 

• ionization: occurs when an incident charged particle deposit energy in a medium 
in such a way that electrons are freed from atoms. It is modeled as a continuous 
energy loss of the incident particle when no secondary emerges from the medium, 
but for sufficient energy ionization results in the creation of 5 -rays, i.e. secondary 
electrons scattered by an incident electron (Moller process) or positron (Bhabha 
process). This is the predominant process by which low energy electrons (below 
lOMeV) lose energy in matter. 

• bremsstrahlung: is the production of photons by a charged particle that is accel- 
erated in the field of atoms. It is the predominant process by which high energy 
electrons (above 10 MeV) lose energy in matter. 

• photoelectric effect: refers to the emission of electrons from the surface of, gener- 
ally, metals, in response to the absorption of an incident photon. This process is the 
main cause of interaction with matter for low energy photons. 

• Compton scattering: a photon can be scattered by a free electron, which gains 
momentum in the process. The scattered photon wavelength differs from the initial 
one by: 

X' — X = — (1 — cos0) 

me 

where 6 is the angle between the incident and outgoing photon. 

The photons produced by the above processes may be involved in further pair produc- 
tion. In summary, the expected processes are all electromagnetic in nature, and the com- 
plexity of the analysis stems from two sources, the stochastic multiple scattering and the 
electromagnetic showering of the electron positron pair. The showering increases with 
the energy of the incident photon. 




PROBABILISTIC EVENT ANALYSIS 


The methodology that we are developing works in principle as follows 

• Run the simulator a large number of times, with the incident photon drawn from 
a uniform distribution. Retain those runs of the simulator where the same detector 
strips fired as in the event. Estimate quantities of interest (means etc. of the photon 
angles) from the runs retained. 

Clearly the system is sufficiently complex that the number of runs is prohibitively large, 
and the description in this section is of a methodology that aims to only simulate those 
runs that will be retained. 

Bayes theorem tells us how to update our beliefs when new data arrive, and how to 
use the distribution from one observation as a prior when analysing the next observation. 
Applied to one event of the LAT, we can impose a time-ordering on the data, and 
incorporate the microstrip firing information one detector layer at a time. The advantage 
of this sequential approach is that the distribution over 6 (azimuth, elevation, energy) 
becomes gradually more peaked as the data is incorporated, making the search for the 
region of parameter space that has significant probability mass more straightforward 
[4,5], 

We adopt ideas from particle filtering, and represent the distribution over 6 by a set 
of uniformly weighted random samples, which are called the “particles ”[6, 7]. Each 
particle represents a possible event - it contains information about the direction and 
energy of the incident photon, and also “nuisance parameters” representing such things 
as where in the detector the pair production occured, the scattering information as the 
electron and positron pass through the instrument, etc. The latter parameters are not of 
primary interest, but are needed to compute the probability of an event. 

Initially we simplify the system to that shown in figure 4. It has been reduced to two- 
dimensions, and we assume that the microstrip detectors are perfect - there are no missed 
detections, resulting in a zero-one likelihood - any simulated event which triggers the 
microstrips that actually fired has likelihood one, and any simulated event that either 
does not fire all the detectors, or fires extra detectors, has likelihood zero. 


Methodology 


The particles at layer i are distributed as 1 and to incorporate the data at 

layer i + 1 requires updating this to p(0 1 . i+1 |<7 1;i - + i), where = {0 1;i ,0 i+1 }, 6 i+1 
being the new nuisance parameters ( e.g . scattering angles) at this level. 


:i+i l^i.-i+i ) - p(d i+l ] 0 f+ i )p{6 l :i \d x :i )p ( 6 i+1 ) 

p(d i+ j | e !+l ) is the likelihood of observing the data at layer i + 1 , and p(0 i+l ) is the prior 
on the nuisance parameters, for example the distribution of the scattering angles defined 


1 We use the subscript 1 : i to denote the set of variables at layers 1 . . . i. 



FIGURE 4. Simulation of the interaction of a 200MeV photon with the LAT. 


by the physics of electrons passing through the material of the detector. 

The standard manner to compute this update is to generate 0 I+1 by sampling from 
p(d i+l ), and then weight by the' likelihood p(6 i+l \d i+l ). However, sampling from 
p( 6 i+ 1 ) result in most of the particles having trajectories that do not intersect the 
detectors that fired, and thus having likelihood zero 2 . Instead, we rewrite the update 
equation as 

P(e l: i+ | K:, +i ) « P W +1 

where p iS (0 i+ j) is an importance sampling distribution [8], that we choose to be con- 
centrated around the values of which result in large likelihood values. The particles 
then have attached a weight given by p{d i+ i[0 i+1 )p(0 i+l )/ p ls (d i+l ) and the particles 
are resampled from the distribution defined by the weights, to give a set of samples that 
represent Mh+i)* 

Note, however, that the elements of 0j. J+1 that are of primary interest, (the pho- 
ton angle and energy), do not change as the data at subsequent layers is incorporated. 
This causes the number of particles that are distinct with respect tO the incident pho- 
ton angles and energy to be rapidly reduced, and the representation of the distribution 
to be poor. For this reason, after each new layer is incorporated, we preform a few 
Metropolis-Hastings updates [9] on all the elements of This leaves the distribu- 

tion unchanged, but reintroduces diversity into the samples. 

Consider the event shown in figure 4. Initially the particles are uniformly distributed 
— the gamma ray may have come from any point in the sky, and entered the instrument 
at any x-position. If we consider the first detector that fires as we move down the stack, 
then it is clear that the pair-production event must have occured in a small region above 
that microstrip, so that the electron and positron pass through the rmcrostnp(s) that fired. 


2 We could instead sample from the likelihood and weight by the prior. While this is possible for the 
simplified detector, this is likely to be infeasible for the actual instrument. 


We therefore construct the first importance sampling distribution to focus the position of 
the pair production event to be in this small region, and the angle of the incident photon 
to be in a small cone determined by the firings at this first level, and those at the next 
level - the size of this cone will depend on the energy of the photon, as that determines 
the distribution of the opening angle. At subsequent levels we construct an importance 
sampling distribution which puts large probability mass on the scattering angles that 
deflect the trajectories such that they intersect with the detectors that fired. 


SIMULATION RESULTS 

We present results from simulations of a simplification of the detector. Simulations in 
2D are shown in figures 5 and 6; a 3D simulation is shown in figure 7. The gamma 
conversion probabilities are those given in [3], as are the distributions of the energy 
split between the electron and the positron (see also figure 3). To simulate the multiple 
scattering, ionisation and bremsstrahlung processes, we used Geant4 [10] to simulate 
the interaction of an electron with a tungsten foil, and fit a gamma distribution to the 
resulting scattering angles for a wide range of electron energies. For the event analysis 
we assumed that the energy split is known - this is necessary in this simplification 
because we do not model the loss of energy by the electron and positron, so there is 
no information about the energy in the scattering angles along the electron and positron 
tracks. 

Figure 5 shows a 200MeV gamma ray normally incident. The energy split is almost 
equal (107:93MeV), and so the “opening angle” is almost symmetric. However, one of 
the tracks happens to undergo more significant scattering than the other. The histograms 
of the photon angle and position as data from subsequent layers are included is also 
shown. The distributions become clearly more peaked. Note, however, that because of 
the finite width of the detector strips, and the fact that the incident photon hit at the very 
left of the strip, the posterior distribution is biased away from the true values - moving 
the interaction point to the right and increasing the incidence angle slightly allows tracks 
that have smaller scattering angles, and are therefore more probable. 

Figure 6 shows another 200MeV photon, where the energy split is much less equal 
(154:46MeV), but where the opening angle is much smaller than the event in figure 5. 
The low energy track undergoes much more scattering (the high energy track in this case 
is almost straight). Because of these characteristics of the event, the parameters of the 
actual event are very close to the modes of the distributions. 

Figure 7 shows a 3D simulation, and the estimates of the distributions of the incidence 
angles and the pair production position. These distributions also become more peaked as 
data at subsequent layers is incorporated. However, after four layers particle starvation 
became a problem, indicating that better importance sampling distributions need to be 
developed. 
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FIGURE 5. Estimated distribution of the incident photon angle and position as data from each layer is 
incorporated (left). The actual trajectory (bold) and a sample of the trajectories represented by the particles 
(right). 
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FIGURE 6. Estimated distribution of the incident photon angle and position as data from each layer is 
incorporated (left). The actual trajectory (bold) and a sample of the trajectories represented by the particles 
(right). 


DISCUSSION 


We have outlined an analysis methodology which we believe to be extremely well suited 
to the analysis of events from silicon microstrip type detectors such as those on the 
GLAST Large Area Telescope. We have shown through (admittedly very simplified) 
simulation the basic applicability of this methodology. Current and future work is de- 
voted to incorporating more of the physics processes into the analysis (especially sec- 
ondary particles, which have the potential to significantly complicate the analysis; and 
energy loss by the particles as they traverse the detector), and improving the likelihood 
model — it is necessary lo incorporate a better model of the firing of the micro strips, 
and include the fact that noise in the electronics may cause firings unrelated to the event 
itself. 
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FIGURE 7. A simulated event for the simplified 3D detector (left). Estimated distributions of the 
azimuth and elevation angles, and the (x,y) position of the pair production event (right). 
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